#include "uhm/common.hxx"
#include "uhm/util.hxx"
#include "uhm/mat.hxx"

#include "uhm/mat/linal.hxx"


using namespace uhm;

int main (int argc, char **argv) {

  initialize(argc, argv);
  
  Mat_Linal_ A(Mat_Base_::DOUBLE, 10, 10);
  A.create_flat_without_buffer();
  A.create_flat_buffer();
  A.create_hier(3,3);

  Mat_Linal_ B(A), C(A), D(A);   
  B.create_flat_without_buffer();
  B.create_flat_buffer();
  B.create_hier(3,3);

  C.create_flat_without_buffer();
  C.create_flat_buffer();
  C.create_hier(3,3);

  D.create_flat_without_buffer();
  D.create_flat_buffer();
  D.create_hier(3,3);

  random_matrix(A);
  random_matrix(B);

  gemm( Mat_Base_::NO_TRANSPOSE, Mat_Base_::NO_TRANSPOSE,
        1, A, B, 1, C );

  std::cout << C;

  linal::gemm( LINAL_NO_TRANSPOSE, LINAL_NO_TRANSPOSE,
               FLA_ONE, ~A.hier(), ~B.hier(), FLA_ONE, ~D.hier() );

  axpy(-1, C, D);

  std::cout << D;

  A.free_all();
  B.free_all();
  C.free_all();
  D.free_all();

  finalize();

  return true;
}


